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We investigate the transport properties of a bilayer exciton condensate that is contacted by 
four superconducting leads. We focus on the equilibrium regime and investigate how the Josephson 
currents induced in the bilayer by phase biases applied to the superconducting electrodes are affected 
by the presence of an exciton condensate in the bulk of the system. As long as the distance between 
the superconducting electrodes is much larger than the exciton coherence length, the Josephson 
current depends only on the difference between the phase biases in the two layers. This result 
holds true in both short- and long-junction limits. We relate it to a new correlated four-particle 
Andreev process which occurs at the superconductor - exciton condensate interface. The system we 
investigate provides an implementation of the supercurrent mirror proposed by Kitaev as a viable 
way to realize topologically protected qubits. 
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1. INTRODUCTION 

Exciton condensates (ECs) are ordered states of mat- 
ter in which macroscopic phase coherence is established 
through the condensation of electron-hole pairs. Since its 
early prediction in the sixties^' there has been consid- 
erable experimental and theoretical activity in this field. 
Signatures of exciton condensation have been reported in 
quantum Hall bilayers^, and in optically-excited excitorP 
and exciton-polaritorP cold gases. 

Early on it was understoocP that spatially separating 
electrons and holes into two distinct and decoupled (semi- 
conducting) layers could be extremely useful to suppress 
electron-hole recombination, and in this way to enhance 
the possibility of realizing an exciton condensate. In the 
quantum Hall regime, exciton-condensation in bilayers 
can be realized in equilibrium by condensing electrons 
and holes in conduction band Landau levels that are lo- 
calized in separate layers^. 

When single-particle inter-layer tunneling processes 
can be neglected, exciton-condensation in bilayers 
is equivalent to spontaneous inter-layer phase coher- 
ence^FI^]. a resu i^ double-layer ECs support a dis- 
sipationless "counterflow" transport channel in which all 
electron-hole pairs drift together, giving rise to equi- 
librium counterpropagating currents in the two layers. 
When the two layers of a bilayer EC are sepa rately con- 
tacted^ remarkable transport anomalie d 3 * 14 ^ -^, associ- 
ated with its neutral counterflow supercurrents^, are ob- 
served. In the absence of a magnetic field, bilayer exciton 
condensation can be enabled by gating the Fermi level in 
one layer to the conduction band and the Fermi level of 
the ot her la yer to the valence band. Two experimental 
g rou p^20j2ij j^yg recently reported the observation of an 



anomalous upturn in the Coulomb drag transresistivity 
as the temperature is lowered. This upturn is interpreted 
as being due to stro ng pa iring fluctuations that precede 
exciton condensation^^ and thus serves as a precursor 
signal for the transition that is similar to the enhance- 
ment of conductivity in superconductors due to super- 
conducting fluctuations above but close to the critical 
temperature^. 

Despite growing experimental evidence, a definite con- 
firmation of exciton condensation in bilayers in the ab- 
sence of a magnetic field is still elusive. The identifi- 
cation of new effects that may highlight the physics of 
ECs is thus highly desirable. In this Article we explore 
hybrid phenomena in which counterflow supcrfluid cur- 
rents are combined with superconducting electrodes. We 
show that when two pairs of superconducting electrodes 
are connected via a bilayer, as sketched in Fig. [T] the 
Josephson current is dramatically influenced by exciton 
condensation. Depending on the bias configurations of 
the superconducting electrodes the device can exhibit a 
variety of behaviors, including the exciton blockade effect, 
where supercurrents are suppressed, and the superdrag 
effect, where the system acts as a perfect current mirror. 

When a superconductor is in contact with a normal 
metal, Cooper pairs can leak across the interface. As 
a result there exists a non-vanishing pair amplitude in 
the normal metal (proximity effect) which implies that 
coherence between electrons and holes in the same layer 
is induced by the coupling to the superconductor. The 
proximity effect is intimately related to the microscopic 
mechanism of transport through superconductor-normal 
metal interfaces. At voltages and temperatures below the 
superconducting gap the dominant process is Andreev re- 
flection^. An electron incoming from the normal metal is 
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FIG. 1: a) Schematic of the four-terminal device studied in 
this work. A double layer is coupled to four superconduct- 
ing electrodes. The system is phase biased and a Joseph- 
son current can flow without dissipation. The bilayer is in 
an exciton condensate phase, b) Assuming translational in- 
variance along the contacts, the system can be considered 
as one-dimensional. In a mean-field approximation the su- 
perconducting and exciton ordered phases can be described 
through space dependent order parameters T(x) and A(a;), 
respectively. The order parameter T is different from zero 
only in a region of length L sandwiched between the elec- 
trodes. To allow for a current-currying state, the phase of the 
exciton condensate must be space-dependent. The four super- 
conducting electrodes have different phases ifi. For simplicity, 
we assume that the amplitude of superconducting order pa- 
rameter is equal to A in all electrodes. 



the emission of a Cooper pair in the bottom layer. 

A brief account of our results was already published 
in Ref. 28 . Here we give details of the derivation of the 
published results and discuss new regimes not analyzed 
previously. The paper is organized as follows. In the next 
Section we consider the case of an ideal interface between 
an exciton condensate and a pair of superconducting elec- 
trodes. We define the model, introduce all the relevant 
energy scales, and discuss the correlated bilayer Andreev 
process which enables supercurrent conversion at the EC- 
superconductor interface. In Sec. |2.3| we discuss the long 
junction limit, providing more details about the results 
presented in Ref. GU Then in Sec. |2.4| we introduce an 
alternative theoretical approach which enables us to ex- 
tend our investigation of the Josephson current to the 
short junction regime. This method is based on the rela- 
tionship between the scattering matrix and the density of 
states. To complete our analysis we study, in Section [3j 
the case of poorly transmitting interfaces. In this case 
the Josephson current can be computed in a perturba- 
tion expansion in the tunneling amplitudes. 

There is a very interesting connection between the sys- 
tem considered in this work and one put forward by Ki- 
taev^ to realize topological protected qubits. In Sec. [4] 
we will explore the similarities and discuss to which ex- 
tent the system considered in this work can implement 
topological quantum computation. In the concluding sec- 
tion we summarize the main results of our work and com- 
ment on the feasibility of testing our findings experimen- 
tally in semiconductor or graphene double layers. Tech- 
nical details of most derivations have been relegated to 
Appendices. 



2. JOSEPHSON CURRENT FOR IDEAL 
SUPERCONDUCTOR-EXCITON CONDENSATE 
INTERFACES 



reflected as a hole at the interface with the superconduc- 
tor, with consequent injection of a Cooper pair into the 
superconductor. Since its discovery in 1966, the study of 
Andreev reflection has offered many surprises. Of partic- 
ular interest is the regime where the (incoming) particle 
and the (reflected) holepreserve their phase coherence 
across the normal metaP^. It is natural to envisage that 
the interplay between phase-coherent electron propaga- 
tion in metals and macroscopic phase coherence in the 
superconductor will be of fundamental importance at an 
EC-superconductor interface, where the neutral EC su- 
perfluid current has to convert onto a (charged) Cooper 
pair current. We will show that this current conversion 
occurs through a new kind of correlated Andreev pro- 
cess. The absorption of a Cooper pair by the exciton 
condensate in the upper layer is always accompanied by 



In the Josephson effeclP^ a supercurrent can flow 
through a weak link in the absence of an electrical bias 
whenever a phase bias is present. In the usual two- 
terminal setup, the weak link can be a tunnel barrier 
(such as a normal metal or a semiconductor sandwiched 
between two superconductors) or any sort of constric- 
tion. In this work we consider^ two coupled weak links 
between four superconducting electrodes, linked by a bi- 
layer of length L, as depicted in Fig. [T] The electron and 
hole densities in the two layers can be tuned separately 
through external gates. We assume that the top (T) 
layer is electron-doped, i.e. negatively charged, and the 
bottom (B) layer is hole-doped, i.e. positively charged. 
The two layers are coupled only through the Coulomb in- 
teraction. (Direct single-electron tunneling between the 
layers is assumed to have been suppressed by a dielectric 
barrier.) 
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FIG. 2: One-dimensional linearized bands for the two lay- 
ers. For a Fermi level location in the figure, top (bottom) 
layer states near the right Fermi point are right-movers (left- 
movers), whereas states near the the left Fermi point are left- 
movers (right-movers) 



2.1. The Hamiltonian 

Assuming translational invariance along the junctions, 
the problem to be considered is one-dimensional (longi- 
tudinal modes with different transverse momentum are 
not coupled). At low energies the dispersion relation 
can be linearized around the Fermi energy e F . In this 
large Fermi-energy limit terms of 0(A/e F ) and 0(T/e F ) 
are automatically set to zero, where A and T are the 
superconducting and EC order parameters, respectively 
(see below) . In the literature the approximation in which 
these terms are neg lected is usually termed "Andreev ap- 
proximation" 1 26 1 31 " ^. 

After linearization the spectrum is split into two dis- 
tinct branches denoted by the two- valued quantum num- 
ber p = ±. On the T (B) layer excitations with p = + 
correspond to right (left) movers and those with p = — 
to left (right) movers (see Fig. [2]). The electron field op- 
erator of layer a — T,B and spin component a can be 
written as 



- ik 



*aa-(aO (1) 



where ^ aa ± (x) are fields related to the Fermi point ±k F a ^ 
that are assumed to be slowly varying over the length 

2tt/4 q) . 

level the exciton and supercon- 
are described by the two order 
{& Ba (x)V Ta (x)) and A(a;) oc 



scales Xp 

At the mean-field 
ducting condensates 
parameters T(x) oc 



(*a^(a : )*at( a; ))i respectively. The EC pairing potential 
pairs holes in the bottom layer and electrons in the top 
layer, while the superconducting order parameter pairs 
right-moving (left-moving) electrons with spin up and 
left-moving (right-moving) electrons with spin down sep- 
arately in each layer. Here we assume that the bilayer 



coherence is between up spins and up spins and equiva- 
lently between down spins and down spins. In principle 
there is a set of equivalent energy bilayer states related 
to these by spin-rotation in one layer only. The mean- 
field Hamiltonian %, quadratic in the fermion operators, 
which captures interactions with these condensates has a 
(generalized) Bogoliubov-de Gennes forrrPll. In Nambu 
notation it reads 



U = 



f + 00 



dx ^2 ^l{x)Hp^ p {x) 
P =± 



where 



H v — 



/—iphv F d x r(x) At(x) 
T*(x) iphv F d x A B (x) 
A* T (x) iphv F d x -T*(x) 



V 



o 



(2) 



A* B (x) —T(x) -iphv F d x J 



(v F being the Fermi velocity) and 



%(x) = (^Ttp(x)^Btp(x),^ 



(3) 



In principle the two order parameters have to be de- 
termined self-consistently. In the regimes considered 
in this paper, however, self-consistency would introduce 
only negligible quantitative changes of the results. A 
more detailed discussion on the situations in which self- 
consistency is unimportant and the order parameters can 
be taken to have a step-like form is discussed in Ref. [35] 
The EC order parameter T(x) is assumed to be uniform 
in amplitude in the EC region (-L/2 < x < L/2) and 
zero otherwise. In order to account for the neutral coun- 
tcrflow current, the phase of the condensate is allowed to 
be space-dependent. Current conservation in one dimen- 
sion for an order parameter with a constant amplitude 
implies the following functional dependence on x: 



T(x) = \T\e li i +q ^ 



(4) 



where the wave vector q and the phase 7 are to be fixed 
by free-energy minimization. The superconducting order 
parameter, on the other hand, is assumed to be different 
from zero only in the electrodes: 

fAe 4 ^, x<-L/2 

0, -L/2 <x<L/2 (5) 

Ae^ T -«, x > L/2 

and 

!Ae^ B - L , x < -L/2 

0, -L/2 < x < L/2 . (6) 

Ae lVB - R , x>L/2 

For the sake of convenience, we have assumed that only 
the phases vary from one superconductor to the other. 
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Because charge is conserved separately in the two lay- 
ers, the number of independent phases is reduced to two, 
which we take to be the phase differences between left 
and right electrodes in the top and bottom layers. 

Because the potentials are constant and are never si- 
multaneously non-zero, it is easy to solve the Bogoliubov- 
de Gennes equations separately in the excitonic and su- 
perconducting regions. (See Appendix|Aj) In addition to 
the distance between the electrodes, L, the other length 
scales which must be considered in coupling different 
regions are the exciton and superconducting coherence 
lengths, £ EC = /w F /|r| and £ s = Hv F /(irA), and, at fi- 
nite temperature, the thermal length L t h = /3hv F [with 
l3=(k B T)- 1 }. 

The Josephson current in this four-terminal system 
has two contributions. The first one is related to di- 
rect quasi-particle tunneling through the barrier caused 
by the gap in the bilayer quasi-particle spectrum. This 
process is present in ordinary Superconductor-Insulator- 
Superconductor two-terminal setups. The second con- 
tribution, which is the main subject of interest here, is 
linked to the possibility of direct dissipationless conver- 
sion of supercurrents into EC counterflow currents and 
can be realized only in a four-terminal setup like the one 
shown in Fig. [T] Importantly, in the limit L ^> £ EC the 
first contribution is exponentially suppressed and only 
the second one survives. In this case, equilibrium prop- 
erties will ultimately depend only on one phase variable. 
We will restrict ourselves to this regime, and address the 
following two distinct limits: 



Long-junction limit L>( s 
Short-junction limit L <C £ s 



(7) 
(8) 



We will further extend our results to the "tunneling 
limit" (i.e. the regime in which the contacts have a low 
transparency) in which the Josephson effect was origi- 
nally studied. Our results in this case are qualitatively 
similar to the ones in the ballistic regime. 



2.2. Four-particle Andreev reflection 

The conversion of supercurrent into counterflow (neu- 
tral) current hidden in the solution of the Bogoliubov- 
de Gennes equations, takes place through a coherent 
Andreev-like process which involves two Cooper pairs in 
the two different layers. We find it useful to give a brief 
illustration of this process before going into the details of 
our results for the Josephson current. A detailed deriva- 
tion will be presented in the following Sections and in the 
Appendices. 

It is sufficient to consider only one superconductor - 
EC interface as shown in Fig. [3j It is useful to imagine a 
thin normal region between the superconducting and EC 
regions in which the bilayer is in its normal phase with 
both order parameters set to zero. We limit our atten- 
tion here to energies smaller than the EC and supercon- 
ducting gaps. As already discussed, the wave function 



has four distinct amplitudes ^ = (ipTf, ipBf, ipTU ^Bi)- 
The superconducting order parameter is responsible for 
the Andreev reflection process in which an electron im- 
pinging on the superconductor from the normal region is 
reflected as a hole and overall current is conserved be- 
cause an electron Cooper pair is created (or equivalcntly 
a hole Cooper pair is destroyed). This process is illus- 
trated on the left side of the bottom layer in Fig. [3j The 
time reversed process is obviously allowed as well and is 
illustrated on the top left of Fig. [3} For the wave function 
this implies the following two relations between up and 
down spin amplitudes in the same layer: 



IpTf — a "0' 



B\ 



(9) 



where a,b are complex number with \a\ — \b\ = 1. Ex- 
plicit expressions for a and b (as well as for the quantities 
c and d defined below) as functions of the energy will be 
calculated later in the paper. An analogous process takes 
place at the normal-EC interface. As shown on the right 
hand side of Fig. [3j a particle (either electron or hole) 
impinging on the EC from the top layer is reflected into 
the bottom layer upon emission (or absorption) of an 
electron- hole pair (exciton). The same process allows a 
particle to go from the bottom to the top layer. The fol- 
lowing two additional relations between top and bottom 
layer amplitudes capture these processes 



(10) 



where c and d are complex numbers with |c| = \d\ = 1. 
Combining the four amplitude relations we conclude that 



"b*d ipT^ 



ac 



b*d=l 



(11) 



Solving Eq. (Ill determines the energy of a state that is 
bound to an interface by the gaps in both superconduct- 
ing and excitonic regions. We refer to this state as to an 
Excitonic Andreev Bound State (EABS in the following) 
It is well known that Andreev bound states in Super- 
conductor - Normal metal - Superconductor systems^J22l 
are responsible for the flow of Josephson current. Al- 
though they differ in origin and in properties, EABS will 
also be crucial for supercurrent conversion from the con- 
tacts to the exciton condensate. Indeed, using general 
arguments based on current conservation, it is natural to 
expect that a phase bias applied between the top or bot- 
tom superconducting electrodes will result in a Josephson 
current equal in magnitude but opposite in sign in the 
two layers. This must be true since only counterflow cur- 
rents can flow deep in the EC. When the total energy of 
the system is minimized to find an equilibrium state, the 
counterflow current generated at the left junction must 
be absorbed at the right. The total Josephson current 
then turns out to depend only on the difference of the 
phase differences across the top and the bottom layers, 
as we show explicitly below. 
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FIG. 3: Four-particle coherent Andreev reflection. - This 
illustration is a schematic for the process responsible for 
the conversion of supercurrent from the two superconduc- 
tors into an excitonic counterflow supercurrent. An artificial 
normal (N) region in which both superconducting and EC 
pairing potentials are zero has been inserted between the su- 
perconductors (S) and the exciton condensate (EC). At the 
superconductor-normal interface an electron impinging from 
the normal region is perfectly reflected as a hole (arrow C) or, 
viceversa, a hole is reflected as an electron (arrow A). Current 
is conserved by emission/absorption of a Cooper pair. This 
is the normal Andreev reflection process. At the right side a 
similar but distinct process takes place. The exciton pairing 
potential has the effect of reflecting electrons and holes from 
one layer to the other (arrows B and D). Since real interlayer 
tunneling is absent the interlayer reflection is always accom- 
panied by the emission of an exciton (electron-hole pair) . Two 
excitons are emitted, one for the process B, and one for D. The 
wavefunctions for the complete path ABCD can be matched 
only at a precise energy which depends on the various phases 
at play and is the energy of the excitonic Andreev bound 
state (EABS) mentioned in the text. The process illustrated 
converts Cooper pairs in the two superconductors into two 
excitons. 



2.3. Long-junction limit 

We start our analysis from the long-junction limit, 
Eq. Q, where the superconducting gap is the largest 
energy scale in the problem. The effect of the supercon- 
ducting leads can be accounted for by means of proper 
boundary conditions. Indeed in the deep subgap energy 
regime e -C A, electrons impinging at the interface with 
the superconductor are totally Andreev reflected with a 
reflection coefficient that is energy independent. The 
boundary conditions at the two interfaces x = and 
x = L for the left and right moving fields [introduced 
in Eq. ([I])] were derived by Maslov et alW* and, for the 
present case, they read 





h(0) = 


-iae iiPc 




(0) , 


(12) 




h(0) = 


-Viae ^ 




(0) , 


(13) 




(L) = 


+iae v 


a\- 


XL), 


(14) 




(L) = 


-iae 1 ^ 


,H<|/t 

at- 


XL). 


(15) 



A convenient way to implement these boundary condi- 
tions is to perform a field folding, namely to express the 



fields ^ aS _ in terms of ^ acr + defined on an extension of 
the system to coordinate — L < x < 0. More specifically 
one defines, for x > 



.(-x) = -iae^ L % L _ (x) , 



at+V ") — 

^ai+(~x) = +iae VaL ^l 



(x) 



(16) 
(17) 



The first two conditions ( 12 ) and ( 13 ) at x = are auto- 



matically satisfied by requiring the continuity of the new 
fields. These definitions in turn imply that 



# QCT+ (-L) = ~e^^-^^ aa+ (L) 



(18) 



For a given layer and spin direction, the two-field model 
(p = ±) can thus be mapped into a one-field model $ a(T+ 
defined on the interval —L<x<L, and satisfying the 
twisted boundary condition (18 1. Consequently the free 



Hamiltonian can be rewritten in the interval —L<x<L 
in the form 



He 



-iahv F & aa+ (x)d x if a(7+ (x) 



Similarly, after the folding transformations ( 12 ) and ( 13 ) 
the excitonic coupling term reads: 



(19) 



-L 



dxT(x) W Ttr+ (x)y Bcr+ (x) +H.c. 



where the folded EC order parameter reads 
f(x)- 



(20) 



i(1+1 x +'Pt.l~<-Pb,l) x <0 



r e - 



The Josephson current through the a = ± = T, B layer 
can be evaluated as 



(I a ) = aev F ^ 

<7=t,4-,p=± 



acrp 



(x)^ aap (x)} 



(22) 



Notice that right-movers are characterized by ap = +1, 
so in the top layer their momentum is located near the 
right Fermi point +/c F , while in the bottom layer is 
located near the left Fermi point — fc F , as shown also 
in Fig. [2j Similarly left-movers are characterized by 
ap = —1. Since all the electrodes are at the same chemi- 
cal potential (no voltage bias is applied), the only contri- 
bution to the current is due to the Josephson term. The 
folded Hamiltonian H = T~Lt + Hb + Hr [see Eqs. (19)- 
p0| ] can be straightforwardly diagonalized in the space 
defined by the boundary condition (18). The supercur- 



rent is then evaluated using Eq. ( 22 ) and can be ex- 
pressed (details are provided in Appendix [B]) as the sum 
of a ground state contribution and a thermal fluctuation 
term 28 : 



la — Ia,GS + la 



TF 



(23) 
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The ground state current reads 



'a,GS 



-a lim 

L y->x 



f e i<*q(x-y) p^(°) 4- ak)e ia< - ki ° ] '+ ak )( x -y) _ e ~ ia( i( x ~y) — a ^ e ia ( k{a) ~ ak )( x ~v) ^ 



while the thermal fluctuation current is given by 

ev F ^ f F(fc(°) - k) 



I a ,TF — ~ 2 



E 

fc<°) 



+ fc) 



1 , e l3(-ahv F q+y/\r\2 + (hv F (fc(i)-fe))2) 1 g /3 (ant;pg+- v /|rp+(fe F (fc(o)+fc))a) 

I 



(24) 
(25) 



In the above equations we have defined 
(0) = (2n + l)vr - _ <^ T + <fB 



2nJ 



2L 



XL 



(26) 



with n and J relative integers. The EC phase-winding 
wave vector is fixed by the condition 



qL = 



where n is a relative integer and 

(p T = if T R - (p T L , 
<Pb = <Pb,r - <Pb,l ■ 

Finally the function F(k) is defined as 



2 + (^ BO )2 



(27) 



(28) 
(29) 



(30) 



i.e. it is a Heaviside-like function, smoothed over a 
length £ EC . We observe that F(k) behaves as 

1 



F(k) 



for k — > — oo 



4(^ EC ) 2 

while for |T| — > tends to the Heaviside step function, 
F(ak) -> 6{ak). 

We recall that the superconducting gap has been set 
as the largest energy scale in the problem and thus it 
does not appear in the expression for the current given 
in Eqs. (24) and (25 1. In the following we will discuss 



the properties of the Josephson current at zero and finite 
temperatures. 



2.3.1. Zero temperature 

In the case of a long junction the function F(k) defined 
in (30) varies smoothly with respect to the discrete fc(°) 
spectrum and one can transform the sum in Eq. (1241) into 



an integral. The counterflow current can be related to q 
which depends solely on the difference of the two phase 
differences because of the condition given in Eq 
The resulting expression, given in Ref. 



(271. 



is 



It/b — ±- 



L V 2tt 



VB 



ip T - ip B e [-tt; vr] . 

(31) 



Eq. (31 1 has the form anticipated previously from general 
arguments. As we will show in Sections |2.4.2 - 2.4.4 it 
is possible to obtain the same result by energy minimiza- 
tion. 

The Josephson current It/b depends on the difference 
between the two phase differences (top and bottom). 
This peculiar phase dependence has several interesting 
physical implications. When the top and bottom junc- 
tions are polarized with the same phase bias (ipr = <PBi 
parallel flow) no supercurrents can flow through the EC. 
In this case the Josephson currents experience an exciton 
blockade. In the opposite case of counterflow phase bias 
(ipT = —^>b) the Josephson current flowing through the 
EC is maximal, with a critical value equal to half of the 
critical current of a ballistic one-channel Superconductor 
- Normal metal - Superconductor (S-N-S) junction. This 
is evidence of the fact that the four-terminal device al- 
lows supercurrent (dissipationless) drag. When current 
flows in one layer due to a phase bias in the same layer, 
a current equal in magnitude but opposite in direction 
flows in the other layer. This is a consequence of the 
perfect conversion of exciton current into supercurrent. 
Eq. (31 1 can then be seen as a perfect drag effect for the 



supercurrent. 



2.3.2. Finite temperature 



The contribution to the current due to thermal fluctu- 



ations is given by Eq. (25). It depends both on the sum 
tpr + ipB and on q [i.e. on the difference <fi — ips, see 
Eq. (27)]. As previously mentioned, under the condition 
£ EC <C L, one can fairly well approximate the sum with 
an integral and the resulting expression is a function of 
q only. Furthermore, under the condition f3\T\ 1 one 
can approximate the Fermi functions in Eq. (25) with 



exponentials obtaining the following expression for the 
total current 



l T/B 



±- 



i-y2^rT sinMg ^ /2) e 



qL th /2 

(32) 

valid in the relevant regime Hv F /L <C k B T <C with 
q fixed by Eq. (27). The first term in square brackets 
is the T = contribution, Eq. (31), while the second 



one encodes the effect of thermal fluctuations and scales 
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as e~^' r '. This means that as long as thermal fluctua- 
tions are dominated by the exciton gap, the ground-state 
current is essentially unaffected by finite temperatures. 
Note that this occurs even when the thermal length L t h 
is smaller than the length L of the junction. This is 
in striking contrast to the case of an S-N-S junction (or 
with the case of two decoupled layers), where the critical 
current is exponentially suppressed^. In the presence 
of the EC, Andreev reflection processes occurring coher- 
ently at the two interfaces transform Cooper pairs into 
the electron- hole pairs of the EC, which are protected 
from thermal decoherence by the excitonic gap. Thus 
in the temperature window hv F /L <C k B T <C |T| the 
EC counterflow channel is responsible for an exponential 
enhancement of the critical current. This effect should 
be readily observable as an anomalous persistence of the 
saw-tooth Josephson current as temperature is increased. 



2.4. From the long-junction to short-junction limit: 
the scattering approach 

The method described in the previous Section, which is 



based on the boundary conditions ( 12 ) - (151 and on the 



folded Hamiltonian, is valid only in the limit in which 
A is the largest energy scale (long-junction limit). We 
now wish to extend our investigation of the current also 
to the regime in which the junction length L is much 
shorter than the superconducting correlation length £ s 
[see Eq. (J7J)]. A different approach is thus necessary. To 
this purpose, we observe that since the Josephson current 
is an equilibrium current it can be computed also from 
the relation 



H dip a 



(33) 



where J-j is the phase dependent term of the free en- 
ergy, and (f a is the phase bias applied to layer a. The 



computation of the current through Eq. ( 33 1 offers the 



advantage that J-j can be evaluated from the knowledge 
of the eigenvalues of the system only (the eigenfunctions 
are not needed^). Namely, one can write 



T3 = — J ale sp(e) , 



(34) 



where p(e) is the density of states of the system and e 
denotes the energy measured from the Fermi energy. The 
calculation of the free energy T3 for the system of two 
Josephson junctions coupled by an EC (Fig. [I]) is one of 
the main results of the present Article. This computation 
enables us to recover the long-junction- limit result pre- 
sented in Sec. |2.3| via an independent method and, most 
importantly, to obtain an expression for the current also 
in the short-junction limit. 

Since the spectrum of the Hamiltonian in Eq. ^ is 
symmetric around zero one can compute the energy den- 
sity between zero and +00. Only the phase-dependent 



part of the density of states is relevant in the calculation 
of the Josephson current. The density of states p de- 
pends on each of the four superconducting phases {<p a ,L 
and Pcr) an d it can be derived in a very simple and el- 
egant way (see Appendix O) fro m the scattering matrix 
S of the junction as follows 
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(35) 



For the sake of comparison, it is useful to recall what 
happens in a S-N-S junction (in the absence of an EC). 
The short-junction regime (where A is the smallest en- 
ergy scale) is by far the simplest to treat, since the 
only phase-dependent feature of the spectrum is a sin- 
gle bound state with energy e < A. In the long-junction 
limitpniM] t ne number of bound states increases linearly 
with L and, moreover, in order to properly evaluate the 
current, one has to take into account the continuum (for 
e > A). The new energy scale |T|, present in our system, 
enriches this picture, because one expects a contribution 
to the free energy due to the counterflow current, which is 
related to features of the spectrum at an energy era |T|. 
Thus, even in the short-junction case, one has to compute 
p{e) at least up to 



2.4-1- Scattering matrix and density of states 

One possible approach to evaluate the scattering ma- 
trix is to generalize the method discussed in Ref. 38 to the 
four terminal case. Here we propose, instead, a different 
approach, based on the following idea. The supercon- 
ducting contacts are not infinitely extended on the left 
and on the right, but truncated to a length M so that 
the superconducting pairing potentials reacP^ 



and 



A B (x) - 



ro, 

Ae ttpT - L , 


Ae ttpT - R , 

u 



0, 

Ae llfiB - L , 
0, 

Ae llfiB - R , 



x < -M - L/2 
-M-L/2 < x < -L/2 
-L/2, <x< L/2 
L/2 < x < L/2 + M 
x > L/2 + M 



x < —M - L/2 
-M-L/2 <x< -L/2 
-L/2 <x< L/2 

x < L/2 + M 



(36) 



(37) 



L/2 < 
x > L/2 



M 



Free-electron plane waves, present in the regions x < 
— M — L/2 and x > L/2 + M, can therefore be used to 
define the total, block diagonal, scattering matrix 



S = 



S+ 
S" 



(38) 
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where S + and S are defined through 



and 



/ ipTf+(d) 
■>pBt+(-d) 

\ ^Bi-{d) 



/ ipTt-(-d) 
ipBt-(d) 

V ipBi+{-d) 



(ipTt+(-d) 
^B\+{d) 
M-(-d) 



^Bi+{d) 



(39) 



(40) 



with d = L/2 + M. Here, if)Ta+(x) and ipB<r-(x) denote 
the wave functions of a free (A = and T = 0) right- 
moving excitation, while ipTa-{x) and tpBa+{x) the left 
moving ones. 

This truncation procedure offers the advantage that 
one can treat the discrete and continuous spectrum on 
the same footing, since now there are available free prop- 
agating channels also with energy < e < A, making 
it possible to define a scattering matrix in any energy 
range. The exact energy density with infinitely extended 
superconducting contacts is found by taking the limit 
M 



+oo in Eq. (|35|) 

truncation approac 



We have also checked that this 



:h reproduces the well-known case of a 
standard ballistic two-terminal S-N-S junction. 

A derivation of the relation ( 35 1 between the determi- 



nant of the scattering matrix and the density-of-states 
can be found in Appendix [C] while in Appendix [D] we 
show that the determinant of the scattering matrix ap- 
pearing in Eq. ( 35 1 can be expressed in terms of the 
transfer matrix T, which is easy to calculate for the 
present case of piecewise-constant potentials. In conclu- 
sion, through Eq. ( p34| ), the free energy can be computed 
numerically. In the following subsections we will focus on 
two relevant limits where analytical results can be found. 



2.4-2. Phase-dependent contribution to the free energy 

We first restrict ourselves to the situation of exciton 
coupling characterized by the inequality |T| 3> A,Et, 
where E^ — hv F /L is the Thouless energy. Since in this 
regime the length L of the bilayer is much greater than 
the EC coherence length, quasiparticle propagation in 
the bilayer is completely suppressed. As a result, the 
density-of-states p can be approximated by the sum of 
two contributions 



p(e,A,|r|) ~p EC (e)+/9 BS (e) 



(41) 



The first term p E c(£) = p(s, 0, |T|) accounts for the EC 
alone, while the second term p BS {e) — p(e, A,oo) is re- 
lated to the superconducting electrodes coupled to an EC 
characterized by an infinite gap. 



Combining Eq. (41 1 with Eq. (34), the free energy can 



be written as the sum of two contributions 



The first term, in the limit of large turns out to be 
(see Appendix [E| 



(43) 



which describes an excitonic supercurrent contribution, 
encoded in the phase winding q of the EC order param- 
eter. 



The second term in Eq. ( 42 ) can be understood as the 



contribution to the free energy due to bound states. In- 
deed, in the limit M — » +oo, p BS (e) consists of two 8- 
functions (see Appendix Jf| , i.e. 



where 



and 



Here 



and 



p BS (e) = 5(e - e L ) + 8(e - e R ) 

qL 



e L = A 



cos 



sr = A 



cos 



2 7 2 



= <PT,L - <PB,L 



PR = VT,R - <PB,R 



(44) 



(45) 



(46) 



(47) 



(48) 



denote the phase differences between left and right elec- 
trodes, respectively. 

The wave functions of the bound states are exponen- 
tially decaying both in the EC and in the superconduct- 
ing contacts. In particular, they are localized at the 
superconductor - EC interfaces and they exist indepen- 
dently one of each other. As a result one obtains the 
total free energy of the system 



Fi{qL,"f,(p L ,LpR) 



2ttL 
A 



(qL) 2 - A 



2 



7 + 



<£r 
2 



7' 



qL 
2 



qL 



(49) 



This is the main result of this Section. 

We now proceed to minimize the total free energy 
J-](qL,j,ipi,,<pii) with respect to 7. We find that the 
optimal value of 7 has to be of the form 



7 



</?L + PR 



(50) 



(42) 



with n a relative integer. This optimal value of 7 phys- 
ically ensures that the currents flowing in the external 
leads have equal value but opposite sign on the left and 
on the right. When ip^ 4- <Pr is altered, only the overall 
phase 7 of the EC responds. 
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Substituting Eq. ( 50 ) in Eq. ( 49 1 we find an expression 



for the free energy as a function of qL only: 



Pj(qL) ee min[^ J (qL, 1 )] = ^r(qL) 2 
7 ZirL 



,'VR-<f>L-1qL\ r 
2Acos I ■ 1 , (51) 



which holds when ipn — tpi, — 2qL £ [— 7r; n]. For values 
of (fR — (pi J — 2qL outside this interval one takes advantage 
of the 27r-periodicity. 



We clearly see from Eq. (51) that the free energy (and 



thus the associated Josephson current) depends only on 
phase difference tpn — <pl, which, by using Eqs. $8 
and ( 48 1-( 47 ) , can also be written as 



(29) 



<PR ~ '-PL — <fT ~ fB ■ 



(52) 



Minimizing Tj{qL) with respect to qL enforces the con- 
dition that the EC counterflow current matches the su- 
percurrents carried by the condensate in the supercon- 
ducting electrodes. In what follows we will drop the tilde 
symbol over J-j for notational simplicity. 

In the remaining parts of this Section we will use 
Eq. (51 1 to derive the Josephson current in the long- and 
short-junction regimes. Before doing that, we wish to 
mention that we have checked numerically the validity of 
Eq. (41) by comparing the three quantities p(e,A, \T\), 
/9 EC (e)7^and p BS (e). The result of this comparison is re- 
ported in Fig. |4j where we plot these three quantities 
versus e. The bound-state contribution to the density-of- 
states Pbs(^), plotted in the lower inset, mainly consists 
of two sharp peaks and is negligible for e > 10 Et (the 
small downward peak is related to the continuum spec- 
trum at energies s > A and disappears in the the limit 
M —> +oo). The EC contribution to the density-of-states 
Pec(£)i shown in the upper inset, is essentially zero up to 
e ~ |T|. It turns out that the density-of-states p(e, A, |r|) 
evaluated numerically from Eq. ( 35 1 , and shown in the 



main panel, coincides with the sum of the data reported 
in the insets. 



2.4-3. Long-junction limit (L 3> £s) 

The current in the long-junction regime has already 
been investigated in Ref. [2"5I with the boundary conditions 

Here we 
based on 



and folding method, as outlined in Sec. |2.3 



provide an alternative derivation of Eq. (31 



the truncation method, introduced at the beginning of 
this Section, and on Eq. (49). This proves the soundness 
of the method. 

Since the long-junction regime is characterized by A 3> 
hv F /L, one can, to a first approximation, neglect the 
quadratic term in the phase winding q in the analyti- 
cal expression for the free energy, i.e. the first term in 
the r.h.s. of Eq. (51). Minimization of the second term 



o 



o 
Q 



600 
500 
400 
300 
200 
100 


-100 



600 




100 105 



110 



200 
100 




1 



i 1 1 r 



j i i i L 



1 2 3 4 5 6 



20 40 60 80 

eL/Tivy 



100 120 



FIG. 4: (Color online) In the main panel we show the density- 
of-states of the system p(e, A, |F|) (in units of 1/Et) as a 
function of ene rgy e (in units of Et) as computed numeri- 
cally from Eq. (35 1 for a random choice of the phases 92^,7, 
and qL and for A = 5 -Et and |F| = 100 Et- Lower inset: 
the contribution due to the presence of the superconducting 
electrodes Pbs(e) as a function of e. Upper inset: the contri- 
bution due to the EC Pec(e) as a function of e. These two 
contributions can be considered as independent in the strong 
coupling limit |T| A,hv F /L. 



in the r.h.s. of Eq. (51 ) with respect to qL yields 



T <Pr - <Pl . 
qL = h nir 



(53) 



with n a relative integer. This is the same result found 
in Eq. (27). There it was derived as a condition for the 
existence of solutions of the Hamiltonian, while here we 
interpret it as a relation enforced by the EABS, which is 
increasingly better satisfied as the ratio A/Et grows. 
Eqs. (50)-(53) can be interpreted as "phase anchoring 



conditions" imposed by the presence of EABS, since the 
phase of the condensate at both interfaces and the super- 
conducting phases are no longer independent. 

The integer n in Eq. ( 53 1 is fixed in such a way to 
minimize the first term of the free energy in Eq. (51). 
The latter thus reads 



irhv F 
2L 



<PR-<PL . ,(<PR-<PL 

- mnt 



2tt 



V 2tt 



where nint(a;) denotes the integer closest to x, or 



2^L 



<Pr - 



(54) 



(55) 



if the phase difference ifR — (p^ is restricted to vary in 
the interval [— 7r; it]. 
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Since an EABS decays exponentially in the bilayer, 
current can flow only through the EC and can be calcu- 
lated as following: 



It/b(q) = ± 



hd(qL) ttL 



a — ^r(qL) 



(56) 



As a result, the Josephson current in a given layer is 
obtained by inserting Eq. ( 53 1 into Eq. ( 56 ) to obtain 



It/b(<Pr ~ <Pl) = ±- 



ev F j tp R - ip L 
2ir 



(57) 



with ipji — ipi, £ [ — 7r; 7r]. Eq. (57 1 can also be derived 



takes the form 



using Eq. (|33|). For the top layer, for example, Eq. (33) 

2e SJj 



h difT 



Using Eq. ( 52 ) one therefore obtains 

ev F ( ip T - ip B 



h 



L \ 2ir 



(58) 



(59) 



which coincides with Eq. (57) 



Note that the current has exactly the same functional 
form of a standard long Josephson junction at T = [see 
Eq. (81 )] albeit with a critical current reduced by a factor 



of two. The significance of the particular combination of 
phases appearing in the argument of the current will be 
discussed below. 



2.4-4- Short-junction limit (L <C £s) 

In the short-junction regime the dominant term in the 
free energy is the EC counterflow current contribution 



(60) 



where e = 2L/^ S -c 1. Minimizing Eq. (51) with respect 
to qL one gets 



T e . f<PR-<PL qL 

QL = - sin — 

y 2 V 4 2 



e . ( <Pr-<Pl \ . 2 \ 
sin I : I + U(e ) 



(61) 



where qL inside the phase argument can be neglected 
since it is a small quantity of order e [this can be easily 
seen from Eq. (61 ) itself]. Eq. (61 1 can be interpreted as 



the equivalent of Eq. (53) in the short-junction limit. 

Substituting Eq. (61 1 in Eq. (51 ) and neglecting terms 
of order e yields 



J"j = -2 A cos 



<PR - <PL 



(62) 



Now, as before, the current can be calculated by cither 



inserting Eq. (61) in Eq. (56) or applying Eq. (33) to 
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FIG. 5: (Color online) The Josephson current as a function of 
Lp T — ipB normalized in two different ways. In the upper panel 
the current is normalized with respect to the critical current, 
while in the lower panel it is expressed in units of ev F /L. 



Eq. ( 62 ) (differentiating with respect to (px for the cur- 
rent on the top layer and with respect to Lpg for the 
current in the bottom layer). In both cases the result is 



Lt/b {}Pt - <Pb 
with tfT — <Pb £ [— 7r; 7r 



, eA . 
±— — sin 

h 



(63) 



2-4-5. Crossover regime (L ~ 

In the intermediate or crossover regime (L ~ ^ s ), the 
Josephson current can be understood as interpolating be- 
tween long- and short-junction limits. In this regime a 
derivation of analytical expressions is particularly diffi- 
cult. A numerical analysis leads to the results summa- 
rized in Fig. [5] It is straightforward to minimize Eq. (51 ) 
with respect to qL, the parameter describing the phase 
winding in the bulk of the EC. As we have seen above in 
Eq. (51 ), after minimization, the final result depends only 
on the phase difference ipu — <Pl = '-Pt — <Pb- This can be 
understood by noting that equal shifts in (fi, and <pr can 
always be reabsorbed in 7, keeping the currents flowing 
in the left and right leads constant. Thus the only rele- 
vant degree-of- freedom left is ipu — p>h and, even in the 
crossover regime, we have the same behavior found pre- 
viously (namely exciton blockade and superdrag) albeit 
with a different current-phase relationship. As shown in 
Fig. [5j the current-phase relation does not change dra- 
matically during the evolution from the short- to the 
long-junction limit. When It is normalized with respect 
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to the critical current (upper panel), the various curves 
differ by a small amount when A changes by three or- 
ders of magnitude (see inset in the upper panel of Fig. [5]). 
Nevertheless, it is important to study the peculiarities of 
the current-phase relation since it can reveal a wealth of 
information about the microscopic processes that influ- 
ence the supercurrenlP2. 

A change in the current-phase relation is also expected 
when the magnitude of the superconducting gaps in the 
four electrodes is different (here we assumed for simplic- 
ity that they are all equal to A). The energy scale in 
en and 6l [Eqs. ( |45| ) and (|46|], and as a consequence 
the supercurrent, will change. Again these will be only 
quantitative effects, the fact that the supercurrent de- 
pends solely on ipn — ip^ is unaltered. 



3. JOSEPHSON CURRENT IN THE 
TUNNELING REGIME 

So far, we have analyzed the behavior of the current 
in various regimes of the junction length, assuming per- 
fectly transparent contacts between the bilayer EC and 
the superconducting electrodes. We will now analyze the 
opposite regime, in which the contacts to the supercon- 
ducting electrodes have low transparencies. In this case 
we employ a tunneling Hamiltonian defined by 



H 



E E 

i=L,Ra=T,B 



with the tunneling coupling given by 



(64) 

K (x')+H.c. 
(65) 

while the superconductors and exciton condensate 
Hamiltonians are expressed in the mean-field approxi- 
mation with order parameters A = |A|e lVi '° and T(x) — 

|p|gj(7+9a;)_ 

The pertubative expansion in powers of the tunneling 
matrix element Ti^ a (x, x') can be set as follows. Let /C = 
/C + Ki where JC = J2i, a %s,i,a + ^ec - fJ>N and K\ = 
Ht- One has to calculate the phase- dependent term in 
the grand-potential 



= -i InZ 



where 



Trie 



-P(H-iiN)i 



(66) 



(67) 



is the partition function and /3 = (fc B T) _1 is the inverse 
temperature. In order to capture all the relevant phase- 
dependent terms one has to expand the free energy up to 
fourth order in the tunneling: 



a = n - 



i 



i(0 2 ) o + i(0 4 ) o 



(68) 



The mean (...) is defined as a mean on the thermal 
ensemble (0) o = Tr(e" 0K °O)/Z O = Tr(p o 0) and 2 
and 04 are imaginary-time ordered integrals 



On 



' rP rP 

/ dr / dT'T T [HT(r)-H T (r')] 
Jo Jo 



(69) 



ctcl 



and 







f P 




04 = 


/ dn.. 


. / dr 4 T t [Ht{ti)...Ht(t 4 )} 






Jo 


Jo 


ctd 



(70) 

By [. . . ] c td w e mean that only fully-connected contrac- 
tions are included. In the following we will focus only on 
the phase-dependent contributions. 

The zeroth-order term £7o gives the energy cost due to 
the phase winding q of the EC order parameter 



2t^L 



(71) 



The contributions coming from the superconducting con- 
tacts necessarily stem from the fourth-order term since 
the superconductors are coupled only through the con- 
densate. The phase-dependent ones involve the anoma- 
lous propagators both for the superconductors and the 
EC: 



1 

24 



(04 , phasc-dopcndcnt / q 



1 

24 



dr\ 



P 



dTAX 



(V^l M ,ci ( T l ) ^t<?2 M ,a 2 ( T 2)> X 
(V'EC,<T 1 ,a 1 (Tl)^ c , CT 3^ 3 (T3))oX 

(V>EC,<x 2 ,a 2 ^2)^0,^4,04 ( T 4)>0 X 

,(74 ,24,CK4 

(r 3 )) + H.c. (72) 



where 



s,t,L,T\ 



\ K = -i",,(ri-T 2 ) }_ \ 1 ip(x 1 -x 2 ) 



*L,T 



IA 



L,T\ 



(73) 



and very similar expressions hold for the other contacts 
(£fc = h 2 k 2 /2m - n). We take the magnitude of the 
parameters |Aj iC( | equal for all the contacts. Only the 
phases vary. For the EC "interlayer propagators" we find 

(V>ec,<7,t(ti, x^ipl b (t 3 , x 3 )) = 



E« 



■ip(x!-X 3 ) 



s 



(74) 



The integrations on the tunneling amplitudes in Eq. ( 72 ) 
are left implicit for simplicity, although their calcula- 
tion can be carried out analytically. Only proper com- 
binations of the indices produce a relevant contribution. 
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Cooper pairing requires o~\ = —o 2 = — 03 = 04, and 
since we are interested in the coupling between differ- 
ent layers ct\ = a 2 = T and a 3 = — B. Finally 
one requires i\ — i 2 = «3 = ii = L, R. Indeed in 
the other case i% — i 2 7^ = one has EC propaga- 
tors connecting contacts on opposite sides of the junc- 
tion, which are exponentially suppressed by the gap |T|, 
though. We thus have two different contributions com- 
ing from terms localized at the two sides. Reshuffling the 
indices is equivalent to a permutation of the integration 
variables, which cancels the 1/4! prefactor. The free spin 
index gives an extra factor two. We take the phase of the 
EC order parameter as constant in the region of the left 
contacts even though it is supposed to wind along the 
junction, i.e. the contact regions are small. This means 
that T = T(0) = |r|e* 7 . Another approximation will be 
to take ti >ot {x, x 1 ) = tS^(x — x') in the contact region. 
This leads to the result 

T L = -^| 4 |A| 2 |lf/(/3, |A|, |r|) cos(2 7 - <p L ) , (75) 



(3' 



where 



EE 



IAI 



(76) 



This integral can be calculated, but its exact value is not 
of importance for the purpose of the present Article. An 
analogous calculation for the right contact produces 



^ = --KHAnr|^(/3,|A|,|r|)x 



cos(2 7 + 2qL - <p R ) . (77) 



Note that in the left-contact contribution [Eq. ( 75 1] the 
value 27 of the EC phase at left side of the junction ap- 
pears, while in the right-contact contribution [Eq. (77l] 
the value 2j+2qL at the right side appears. These phases 
are coherently coupled with the difference of the left and 
right superconducting phases ipi, and tpn, respectively. 
Together with the EC bulk contribution this produces 
a free energy, which is in form very similar to the one 
studied in the case of ideal interfaces, modulo a differ- 
ent functional dependence on phase [cos(<p), typical of a 
tunnel junction, instead of cos(y>/2)]. As we said earlier, 
the value of the prefactor is not so relevant here since the 
main result that can be derived from the total free-energy 
functional, 



T = F EC (qL) + Jl[2 7 - (<p L>T - <p LiB )] 
+J r B\2-f + 2qL - (<p RtT ~ <Pr,b)] , 
is that the Josephson current is of the form 



I 



T/B 



= ± I r sin 



/ (fir - <fiB 

V 2 



(78) 



(79) 



for ifT~<PB £ [—Ti - ; 7r] and then extended by periodicity. 
The derivation is very similar to the one given for the 
ballistic short junction (Sec. 2.4.41. 



4. TOPOLOGICALLY PROTECTED QUBITS 

A very important consequence of the fact that the 
Josephson energy depends only on the difference ipx — <£b 



[see Eqs. (31) and (63 1 for long and short junctions, re- 
spectively is the appearance of an almost exact dou- 
ble periodicity in the energy of the circuit illustrated in 
Fig. [6j The double periodicity is analogous to the one 
suggested originally by Kitaev22 for Josephson junction 
arrays with a similar property, i.e. energy dependence 
on the variable ipr — <£>b only. The cross connections in 
Fig. [6] ensure that fluctuations which change 4>t~ 4>b are 
suppressed, leaving 4> — ipx — —<Pb as the only degree of 
freedom. An energy that is function of the phase differ- 
ence ipT ~ ^b may be expressed as F((pr — ^b) — F{24>) 
which is doubly periodic compared to a usual Josephson 
energy. 




FIG. 6: A four terminal device studied which realizes a 
protected qubit. The connections between pairs of ter- 
minals A, A' and B,B' ensure that the phase differences 
(fx ~ <fi(B') — <p(A) and ifB = <fi(A') — f(B) are opposite: 
ifiB = — <£t- This results in a doubly periodic free energy 
dependence on the phase difference </> — <p(B) — <p(A) and for- 
mation of two degenerate states that can be used for quantum 
computation. 

Exact double periodicity of the free energy is 
quite generally promising for decoherence-free quantum 
computation 42 . It implies that the energy of the resulting 
circuit has two distinguishable minima located at <f> = 
and (f> = tt that are separated by a maximum at <f> = ±n/2 
with energy E 2 . The two quantum states |0) and \w), cor- 
responding to phase differences = and = 7r, can 
be used as the logical states of a quantum computation. 
The fact that their energies are equal implies absence of 
the dephasing processes. The large barrier between them 
implies that there is no tunneling and thus no decay. 

Decoherence-free quantum computation becomes pos- 
sible only if one can satisfy a number of constraints on 
the energy scales that we discuss below. We now show 
that these constraints might be easier to satisfy for the 
four-terminal device discussed in this Article than for the 
original "Josephson-current mirror" suggested by Kitaev. 

The height of the barrier is E 2 = irhv F /(2L) for a long 
one-dimensional contact. In general, for a more realis- 
tic two-dimensional system, e.g. a graphene sheet with 
transverse size d, we expect E 2 — r y(k F d)hvF/ L where 
7 ~ 1 is a numerical coefficient. As we explain be- 
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low, the optimal value of E 2 is in the range 1 — 10 K, 
which is achieved for one-dimensional contacts of length 
L = 1 - 10 fim. 

The energy difference between the two minima of 
the potential is due to single-electron tunneling through 
the contact which has been neglected up to now in 
the discussion because it is exponentially small: 5E ~ 
E 2 exp(— L/£ EC ). Because E 2 decreases slowly with the 
increase of the length of the contact and is expected 
to be proportional to the width of the contact for two- 
dimensional structures it should be possible to produce 
a circuit characterized by energy E2/k B > 1 K and 
very small splitting SE/k B < 10~ 10 K. Almost ex- 
act double periodicity of the energy of the circuit im- 
plies that the two quantum states, |0) and \n), remain 
coherent for a time that is limited by the shortest of 
the inverse splitting time h/SE > 1 s and the time re- 
quired for quantum or thermal tunneling across the bar- 
rier. The quantum tunneling is due to the charging 
energy, Eq = e 2 /2C, of the mesoscopic device. This 
energy can be made small by attaching a capacitor as 
shown in Fig. [6] which makes the tunneling amplitude 

t ~ £j /4 #c /4 ex P [ - (n/8)y/E 2 /Ec] of the same or- 
der as the splitting between the maxima. For E 2 > IK 
the capacitor should be of the order of C ~ 1 pF. 
Thermal tunneling across the barrier is exponentially 
small at the typical base temperature of T ~ 20 mK: 
Tth ~ Up 1 exp(— E 2 /k B T) ~ 1 s where io v is character- 
istic frequency of the quantum oscillations within each 
minima that is due to charging energy of the device: 

h^ v = (±l%)^[E2Ec- 

Thus, a device of this type should preserve coherence 
between the two quantum states for a time of the order 
of 1 s. Two quantum states can be controlled by includ- 
ing additional SQUID loops in the connecting loops that 
vary the energy of the logical states depending on mag- 
netic field, as discussed by Kitaev^. The operations on 
the quantum states do not excite quasiparticles in the su- 
perconductors provided that the energies remain smaller 
than the superconducting gap. This condition translates 
into the requirement that E 2 < A. On the other hand 
one wants to keep E 2 as large as possible to suppress 
thermal excitations. These two conditions together im- 
ply that the range E 2 /ks ~ 1 — 10 K mentioned above 
is optimal. 



5. CONCLUSIONS 

In this Article we have calculated the Josephson cur- 
rent between two pairs of superconducting terminals cou- 
pled by a bilayer electron system that hosts an EC. We 
have focused on the regime of strong exciton coupling 
where the bilayer gap |F| is the largest energy scale. 
In this limit quasiparticles cannot propagate through 
the bilayer and the Josephson current is entirely due 
to the conversion of Cooper pair current into counter- 
flow excitonic supercurrent. We have considered both 



the short- (L <C £ s ) and (L 3> £ s ) long-junction regimes. 
The Josephson current, at zero temperature, is given by 
Eq. (63) and Eq. (|3l|) in the short- and long-junction 



regime, respectively. Results for the intermediate regime 
are plotted in Fig. [5] The Josephson current, calculated 
numerically, is plotted in Fig.[5]as a function of the phase 
difference ipx — <Pb for different values of the ratio A/E^ 
spanning both regimes. 



Surprisingly, these two results [Eqs. (63) and (31)] can 



be simply obtained starting from the standard result for 
the Josephson current I(<p) for the two separate (top and 
bottom layer) junctions and applying the substitution 
ip — > (lpt — <Pb)/2- Namely, for a standard S-N-S ballistic 
short junction at T = 



iA 



sm 



(1) 



whereas for a S-N-S long junction at T = 



(80) 



(81) 



In particular, in the short-(long-) junction case the criti- 
cal current is reduced by a factor l/y/2 (1/2). 

The reason why the superconducting phases enter only 
the Josephson current expression in the combination 
{if T — ip B )/2 may be traced back to the relation between 
the superconducting contact phases and the EC phase, 
which is energetically forced by the EABS [as an example, 
see the "phase anchoring conditions" in Eqs. ( 50 )-( [53]) 
or (50)-(61)]. From the derivation given in Appendix IFI 



it emerges that EABS encode a process of correlated 
Andreev reflection, which is the only current conversion 
mechanism available in the presence of a strong exciton 
condensate. For each Cooper pair absorbed in the top 
layer by the condensate a Cooper pair is emitted on the 
same side in the bottom layer. This is the only way for 
a Cooper pair to enter the condensate. 

The factor two in the denominator of the expression 
(ip T — (Pb)/2 can be explained as following. We begin 
by observing that a 27r-change in the phase difference 
corresponds to the tunneling process of a Cooper pair 
from one side of an ordinary Josephson junction to the 
other. The number of Cooper pairs is not conserved and 
the state of the system is unchanged, i. e. it has the same 
free energy and current. For our four-terminal junction 
this process is forbidden since the motion of a Cooper pair 
in one layer must be balanced by a counterpropagating 
pair in the bottom layer. Indeed, if we change ipr to 
tpx + 27r we put the junction in a different current state. 
However, a simultaneous change of <pg by 2tt leaves the 
junction state unaltered. Thus the fact that the current 
is a 27r-periodic function of (ipx — </?b)/2 and not of e.g. 

~ <Pb is a consequence of the fact that electrons are 
transferred through such a hybrid junction in groups of 
four. 

The peculiar phase dependence we have found has sev- 
eral interesting physical implications. When the top and 
bottom junctions are polarized with the same phase bias 
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ipx = <Pb (parallel flow) no supercurrent can flow through 
the EC. In this case the Josephson currents experience 
an exciton blockade. In the opposite case of countcrflow 
phase bias (ip? = ^^Pb) the Josephson current flowing 
through the EC is maximal, with a critical value equal 
to half the critical current of a ballistic one-channel S-N- 
S junction. The four-terminal device allows to realize a 
superdrag effect, i.e. a drag of dissipationless currents. 
When current flows in one layer due to a phase bias ap- 
plied to that layer, a current equal in magnitude but 
opposite in direction flows in the other layer. This is a 
consequence of perfect conversion of exciton current into 
supercurrent. 

At finite temperatures, as long as the EC gap is larger 
than k B T the current is essentially unaffected by ther- 
mal fluctuations. Note that this occurs even when the 
thermal length L t ^ is smaller than the length L of the 
junction. This is in striking contrast with the case of an 
ordinary S-N-S junction (or with the case of two decou- 
pled layers), where the critical current is exponentially 
suppressed^ when L t h <C L, due to thermal decoherence 
affecting a single Andreev-reflection process. In the pres- 
ence of an EC, Andreev reflection processes coherently 
occurring at the two interfaces transform Cooper pairs 
into electron- hole pairs of the EC, which are protected 
from thermal decoherence by the excitonic gap. Thus 
in the temperature window hv F /L <C k B T <C |T| the 
EC counterflow channel is responsible for an exponential 
enhancement of the critical current. This effect should 
be readily observable as an anomalous persistence of the 
saw-tooth Josephson current as temperature is increased. 



5.1. Realization of the device 

Let us now discuss about possible implementations of 
the proposed setup. Electron-hole bilayers have been re- 
alized in semiconductor (e.g. GaAs) double quantum 
wells separated by a thin (e.g. AlGaAs) barrier. In 
these systems it is possible to selectively contact one 
layer by depleting the other one through suitable gat- 
ingSU This technique has been successfully applied in 
"Coulomb drag" experiments, which have recently pro- 
vided indications of EC formatiorp2H22l. I n these systems, 
however, it may not be easy to contact superconducting 
electrodes and to measure equilibrium currents, because 
a large normal gap would arise at the depleted layer be- 
tween the superconductor and the EC. So far, the obser- 
vation of exciton condensation under equilibrium condi- 
tions has been achieved only in quantum Hall bilayers at 
total filling factor = 1- However, quantum Hall sys- 
tems necessarily have current-carrying gapless channels 
localized at their edges, which may alter the physics we 
have discussed above. 

Another possibility is to realize ECs by employing two 
closely-spaced decoupled graphene layers^, hosting a gas 
of spatially separated electron-hole pairs. As compared 
to semiconductor bilayers, such ECs have been predicted 



to exhibit extremely high critical temper at ures^HUl and 
much larger electron and hole densities. Moreover, the 
small distance between the carrier layers, the weaker 
dielectric screening, and the linearly-dispersive bands 
help to increase both interaction and disorder energy 
scales. Furthermore, graphene bands are nearly per- 
fectly particle-hole symmetric, guaranteeing the nesting 
between the Fermi surfaces of the conduction and va- 
lence bands, and favoring the emergence of a coherent EC 
state. Josephson currents flowing through graphene con- 
tacted to superconducting electrodes have already been 
observed by several groups^. For these reasons, the 
use of double-layer graphene seems to be an extremely 
promising direction. Due to their ultra-high mobili- 
ties (even at room temperature), double-layer graphene 
sheets embedded in a matrix of BN layers®! are at the 
moment the most promising graphene-based candidates 
for the observation of EC in the absence of an external 
magnetic field. Inter-layer tunneling can be suppressed 
by interposing a sufficient number of BN layers between 
the two graphene sheets. 

A third realization scheme could be based on 3D Topo- 
logical Insulators (Tlsj^l These recently discovered ma- 
terials exhibit Dirac-like conducting surface states sepa- 
rated by an insulating bulk. Recently, it has been put 
forward that ECs could be realized, at least in princi- 
ple, by oppositely gating the surfaces of a TI thin filnPH, 
or by inserting a thin insulating layer between the top 
surface of a TI and the bottom surface of another Tp21. 
The problem of contacting the surfaces to superconduct- 
ing electrodes has not been addressed yet. Nevertheless, 
this type of implementation may become realistic in the 
near future, in view of the rapid technological advances 
stimulated by the topological protection offered by these 
materials. 



5.2. Possible applications 

The unique properties of the conversion of exciton onto 
Cooper pair supercurrents can be exploited for a num- 
ber of possible applications. In Section [4] we discussed 
how the appearance of an exact double periodicity of the 
free energy in the circuit shown in Fig. [6j allows one to 
realize topologically protected qubits. It is also possi- 
ble to imagine a device in which two superconducting 
electrodes contacted to (say) the top layer are enclosed 
to form a ring-shaped rf SQUID geometry, so that the 
phase difference ipx is directly related to the magnetic 
flux tpT — 27r$/$o + 27rn. In response to a magnetic field, 
an induced Josephson current It flows in the top layer 

-It 



and, according to Eq.(31 1, an opposite current Ib 



flows in the bottom layer. Whenever the magnetic field 
changes the flux by a fluxon <l>o, the currents in both lay- 
ers are reversed. Current sign switches detected in the 
bottom layer count the fluxons present in the top layer 
ring. If the magnetic flux is generated by a monotonous 
analog input signal, the system effectively converts it into 
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a sum of current switch pulses, i.e. to a digital signal. 
The system is therefore an analog-to-digital converter. 
A generalization to non-monotonous input signals can 
easily be achieved by using two double junctions. An- 
other possible application is photodetection. Indeed the 
excitonic current is altered if the layers are exposed to 
an electromagnetic or noise source, resulting in a modi- 
fication of the Josephson currents in the two junctions. 
Finally, we also observe that if electrodes in both lay- 
ers are enclosed in a ring-shaped geometry, any excitonic 
supercurrent reversal generates a fluxon-antifluxon pair, 
indicating that this device can realize Josephson fluxon- 
antifluxon transistors. 
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Appendix A: Solutions in uniform pairing potentials 

In order to calculate the energy density one needs to 
calculate the transfer matrix which in turn can be con- 
structed from the eigenfunctions obtained by piecewise 
solving Hamiltonian ([2| in the separate regions where 
the pairing potentials are constant. 



1. Exciton Hamiltonian 

Let us start with the excitonic Hamiltonian in the mid- 
dle region 

(—iphv F d x T(x) \ 
T(x) = \T\e^ + ^ 
T*(x) iphv F d x J 

(Al) 

Consider first solutions with energy \s — phv F q/2\ > \T\ 
and wavevector k. Introduce the following very useful 
notation 



cosh 0=i/l 



sinh 9 = 



hv F k 



hv F k 



(A2) 
(A3) 



re -pr0/2 e -i(-y+q X )/2 i 



Then the solutions can be written 

(l>£\xj\ / eVrBI1 e i{n+ q x)l1 

where r is a sign defined as 

r = sign (e — phv F q/2) 
They have dispersion 



e(k) = phv F q/2 ± v/|r| 2 + {hv F kf . 



(A4) 
(A5) 

(A6) 



Let us now analyze the case \e — phv F q/2\ < \T\. These 
solutions can be obtained from the ones found previously 
upon substituting 



k — > in and 



and updating the definitions ( A2 1 and ( A3 ) 



cos = 1/1 — 



hv F 



Hv f k 



(A7) 

(A8) 
(A9) 



We should also include the corresponding solutions for 
the down-spin wavefunctions, but they can be obtained 
simply by noting that 



H 



BG,4. 



-n 



EC,t ' 



(A10) 



2. Superconducting Hamiltonian 



The Hamiltonian that describes a bulk superconduct- 
ing contact in the top layer is 



-iphv F d x Ae tv 
Ae~ lv iphv F d x 



(AH) 



The solutions are almost identical to the ones written 
previously, so we are just going to fix the notation. We 
need only the solutions with positive energy. Using 



cosh (8 = y 1 
sinh f3 e 

we have the solutions (e > A) 



hv F k s 
A 

Hv F k s 
A 



(A12) 
(A13) 



e i V /2 e p/3/2 



1pTp(x) 



e -iip/2 e -pf}/2 
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The dispersion is 



:(k s ) = VA 2 + (hv F k s y 



In the case < e < A the evanescent wave solutions 
are again obtained with the substitutions 



k a — > in s and j3 — > i/3 . 
from which it follows 



1 



sin/3 



Hv f k s 
A 

Hv F K s 



(A16) 

(A17) 
(A18) 



When considering the bottom layer superconducting con- 
tacts one needs to note that the Hamiltonian is given by 
-P- 



(A15) Eq. (All) withp 



Appendix B: Diagonalization of the Folded 
Hamiltonian 



We now give a derivation of equation ([25]). Start with 
the field operator written in the form (using the folded 
Hamiltonian formalism) 



(k+q/2)x~ 



C-kr 



l 



r=± k 
\m-f -i(tp T ,L-<PB,L)/2 



J(k-q/2)x; 



r=± k 



(Bl) 
(B2) 



m-y is defined by (^t,l ~ <^b,l)/2 — 7 = m 7 7r and the sign r = + corresponds to the upper branch (positive energy 
with respect to e P ) of the dispersion, while r = — to the lower branch. Under the condition hv F q <C |T| one has that 
the occupancy is purely determined by the sign r of the branch, namely (cl k c r k) = 5 Tt - and (cVfccL.) = <$ r ,+- The 
expression for the current in layer a = ± reads 



I a {x) = aev F (*L+(*)*a<7+(a:) - *L-(*)*a<r-(a 
ff =t4 



(B3) 



Using the folded Hamiltonian one can rewrite 



I a (x)=aeV F (*L+(i")*aa+(3;) - *aa+(-aj)*L+(" 



<r="U 



5^ lim (*L+(^)* QCT +(y)-*a CT +(-y) , J'L+(- 



(B4) 



where in the last line we have introduced the correct point- splitted definition of the product of two field operators. The 
evaluation of the (free) energy and consequently of the supercurrent has to be performed by adopting a point-splitting 
procedure in order to resolve the ill-defined product of the field operators at the same point in space. This requirement 
is necessary because of the linearization procedure around the Fermi points. Introducing now the solution for the 
fields (Bl) and (B2) one has 



Ia(x) = (I a (x)) = aev F y lim 

' ^ 11 — iff 



e i{k+aq/2){x-y) c kr ) - e ii - k+aq ^ 2 ^ x ^ y) (c kr c\ 



r=± k 

a— lim V F(-ak) f e -^+«e/2)(x-t/) f, E - ^ _ e i{k+*q/2){x-y) ff_E-(k))) 
F(ak) (e- l{k+aq/2){x - y) f{E+{k)) - e l{k+aq/2){x - y) f(-E+(k))^ 



(B5) 
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where / is the Fermi function and E r (k) — hv F q/2 + r^J (hv F k) 2 + \T\ 2 the dispersion. Manipulation of this latter 
expression allows to single out the factors related to T = and finite temperature fluctuations 



I a (x) = —a 



lim V \F{ak)e l(k+aq/2)i - x - v) - F(-ak)e- ilk+aq ^ x -^ 

k 

+2a e -^ J2 [F(ak)f(E+(k)) - F(-ak)f(-E-(k))] 



(B6) 



In the folded geometry the allowed wavevectors have the following form (J is an integer) 



k = fe (0) + k jfc(°> = 

2X 



(fT + <Pb + 2nJ 
4L 



Let us study first the ground state contribution 



(B7) 



I a ,Gs(x) = -a— lim 



^ e iaq(x-y)/2 p^iO) _|_ ^ e ia(k (0) +uk)(x-y) _ & — iaq(x-y)/2 p^(0) _ fc) e ia ( fe( ° 3 -ak){x-y)\ 



Lfe(o) 



(B8) 

in the case that ^g^, 1 ^> Afe = Ao = 2n/2L = tt/L, or in other words ^ ex <C L, the function F(/c) appearing in the 
above equation varies smoothly with respect to the discrete fc*- -* spectrum and one can transform the sum into an 
integral 



ev ¥ 



+00 



dfc (0) (e la( i( x -y)/ 2 F(k {a) + k^ e Mk m +ak)(x~y) _ e -iaq{x-y)/2p^ k (Q) _ ^gi«(& (0) -ak)(x-y)\ 



L y^x 

(B9) 

Introducing new integration variables p — k^ ' + ak in the first term and p = k^ — ak in the second term, we notice 
that the dependence on k (i.e. on the sum tpj* + ipg) disappears, whereas the dependence on q (i.e. on the difference 
ip T — ip B ) remains. One thus obtains 



I a ,Gs(x) = -a— lim 

7T y— >x 

= —la lim 

7r y— >x 



Jaq(x-y)/2 



sin[q(x - y)/2] 



(x - y) 



d P F(p)e iap{x - y) 
dF(p) 



-iaq{x~y)/2 



dp F(p)e 



iap(x-y) 



dp - 



Op 



f e iap(x-y)\ 



(BIO) 



We notice that since F(p) is bounded < F(p) < 1, the integral of its derivative converges, and to lowest order in 
x — y one can set 



dp dF (P) e i<xp(x-y) 
dp 



Thus one obtains 



ev 

Ia,as(x) = -2a— lim 

7T y^fX 



+ °° dF(p) 

dp - F(+oo) - F(-oo) = 1 

-co OP 

sva[q(x-y)/2} 

{x - y) 



(-1) 



= a q , 

7T 



(Bll) 



(B12) 



which (as expected) does not depend on x. Recalling that q is given by (27) and the integer n must be chosen in 
order to minimize the free energy oc q 2 , one has the current in the limit L 3> £ E c &s given in Eq. (31). 
Consider now the finite temperature contribution 



I a ,TF(x) = -2"^- X! x 

fe(0) 



F(fc(°) - k) 



+ k) 



1 + cxp 



P (-ahv F q/2 + y/\T\ 2 + [hv F (k(°) - k)}- 



1 + exp 



(i (ahv F q/2 + y/\T\ 2 + [hv F (k(°) + k)}- 



(B13) 



We notice that the above expression depends both on k (on the sum Lp? + ifg) and on q (the difference ifT — Vs)- 
Under the conditions |F| ^> hv F q/2 we observe that the sign of the exponents appearing in the Fermi functions is 
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fixed and independent of k^°\ The means that (regardless of the temperature) the integrand varies smoothly with 
k^°\ and one can fairly well approximate the sum with an integral, obtaining 



I a .TF{x) = -2— / ( 

n J-oo 



F(fc(°) - k) 



F(fc(°) + k) 



1 + exp 



p (-ahv F q/2 + y/\T\ 2 + [hv F (kW - k)] 2 } ] 1 + exp [/3 (ahv F q/2 + ^\T\ 2 + [hv F (k<-°) + fc)] 2 ) 



(B14) 



Each of the two terms separately converges due to the finite temperature. One can than introduce new integration 
variables p = fc(°) — k in the first term and p = kS a ^ + k in the second term. In doing that, we thus obtain 



,TF 



(x) 



+00 



dp x 

oo 

F(p) 



F(p) 



1 + exp 



/3 (-ahv F q/2 + ^|r| 2 + (^ F p) 2 )] 1 + exp U (ahv F q/2 + y/\T\ 2 + (hv F p) 2 ^j 



Furthermore under the condition /3\T\ 3> 1 (£ E c *C L th ) we can approximate 

1 -f}{±ahv F q/2+^\T^ + (hv F p)^ 



1 + exp 



/3 ( ±ahv, 



, q /2 + V|r| 2 + (hv F p)^ 



obtaining 

Ia,TF(x) 



,ev F 



d P F(p) 



-p(-ahv F q/2+y/\r\* + (hv FP )^ _ -0 (a hv F q/2+y/\r\* +(hv F p)^ 



= -4— sinh[pahv F q/2} f °° dpF(p) e -^\ r \ 2 +( h ^P) 2 

K J-oo 

Recalling the definition of the function F(p) 



F(p) = 



i vlrp + (hv F p) 2 ) + hv F p 



2 v/|r| 2 + (hv FP ) 2 ) 

and noticing that the integrand is even in p, one can write 



t r \ A eVF ■ ur« h /oi A ^ a/ r 2 + (fe F p) 2 ) + ?w F p _ /J ./|r|2 + ( toFP) 2 
The integral can be carried out in terms of the Bessel function Ki 



eiri 

I a TF(x) = -4a- 1 — - smh(0Hv F q/2) 
irn 



EC - i^t- 

- — e ? ec 



To first order in £ EC /L th and using the asymptotic espansion of the Bessel function 

Ki 



£ec 



e £ ec 



2L t 



(B15) 



(B16) 



(B17) 



(B18) 



(B19) 



(B20) 



(B21) 



one retrives Eq. ( 32 1 (where the ground state contribution has been added) . We stress that this result holds as long 
as £ec <C L 1 L th . ^ 

Appendix C: Prom the scattering matrix to the two ingoing channels (01,02) and two outgoing channels 

energy density 



The aim of this section is to offer a simple derivation of 
Eq. ( 35 ) . Consider for simplicity a simple scatterer with 
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(bi, b 2 ) described by a 2 x 2 scattering matrix S 



Appendix D: Scattering matrix determinant 



ri t 2 \ ( ai 
h r 2 ) \a 2 



(CI) 



The coefficients j-j and t% are arbitrary function of energy, 
constrained only by unitarity. Now suppose that on both 
sides of the scatterer there are perfectly reflecting mirrors 
enforcing the following conditions on the wavefunction 
amplitudes 



b 2 = a 2 e l 



(C2) 
(C3) 



cj) is an arbitrary phase that can be changed at will. Thus 
the allowed values of the energy can be found solving the 
following equation in e 



i(f>l{e) p i(f>2{e) 



J<t>\ — . 



9 ) 



(C4) 



are the eigenvalues of the scattering matrix. 
Eq. ( C4 ) has in general many solutions since it is satisfied 
when 



4>i(e) = ct> + 2kir, £=1,2. 



(C5) 



Take (f> + 2kn of the form 2irn/N with neZ and with 
N sufficiently big to have a fine sampling of <fii(e). The 
number of states in the small interval [e; e + Ae ] is given 
by 



1 9 ^ A AT 

2tt ae 



(C6) 



This expression has a well defined limit for Ae — > and 
N — > +00, with their product kept costant, and is pro- 
portional to the density of states. In order to determine 
the proportionality factor, one can consider the case of a 
trivial scatterer, i.e. a region of free propagation of a cer- 
tain length. The right energy density is given by fixing 
AeN = 1. In the end one has 



P(e) = 



1 (d<j>x{e) , dfa(e) 



2tt 



de 



(C7) 



an alternative 
Eq. pi. 



way of writing the more general formula 



We first calculate the block diagonal transfer matrix T 
of the system 



T 



T+ 
T- 



(Dl) 



where T + and T are defined by the expression (d 
L/2 + M) 



ipTip(-d) 
V ipBip(-d) 



(D2) 



ipBtp(d) ^ T p 
ipTip(d) I 
V tpB\p{d) J 

where ip aap (±d) are the channels at the two ends (±rf) 
of the device (with superconducting contacts truncated) . 
By dividing the system into three regions, one can write 
that TP = TlxTgxTl^, where T S P R is related to the right 
superconducting leads, T^ L is related to the left super- 
conducting leads, and T£ x is related to the exciton double 
layer. 

Once the (unnormalized) eigenfunctions for the three 
regions are given, their transfer matrices can be written 
in the formal 

/ e ik(x a -xi) \ 

22x2(^2,^1) = A(x 2 ) I Q e - ik ( X2 -x 1 ) J A 

(D3) 

(2^(l,r) is block diagonal in the layer space and similarly 
Te X is block diagonal in spin space, the above formula 
is relevant for each 2x2 block. x 2 ,xi are the two ends 
connected by the transfer matrix). 

A is the matrix of the wavefunction components re- 
ported in Appendix [A] Considering for instance the up- 
per block of Tlx (acting on V'Tfi V'Bt s P ace ) one has 



A(aO 



oP re/ 2 i(-f+qx)/ 2 -pre/ 2 i(~f+qx)/ 2 



re -pre/2 e -i(-y+qx)/2 re pr9 /2 e -i(-i+qx) /2 



(D4) 

In particular, for the left (right) superconducting con- 
tacts region of length M we find (for the full 4x4 ma- 
trices) 



SMR) S inh (3 



( sinh (J3 + ipMk s ) 

wah{fi- ipMk s ) 

ipe-WT.nm sin (Mk s ) 



-ipe itpT MR) sin (Mk s ) \ 

ip^fBMR) sin (Mk a ) 

sinh(/3 - ipMk s ) 



V 







-ipeT %i CB,L{R) s in(Mfc s ) 







sinh (/3 + ipMk s ) J 

(D5) 
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while for the middle region of the EC double layer of length L we find (6 + , k + is calculated with p = +) 

T p = x 
EC sinh0 p 

/ e 4 « L / 2 sinh (0 p + iprfcpL) -ipe** sin (fc p L) 

ipe" 47 sin (fcpi) e~ l i L / 2 sinh (6 p -iprk p L) 

e -wi/2 S i n h (6» p - iprkpL) -ipe^ sin (fc p X) 

ipe^ sin (k p L) e iqL/2 sinh (6 p + iprk p L) J 



(D6) 



where # p and fc p are different functions of the energy in 
each branch 



cosh0 p = 



sinh#p = 




k„ = - — sinh 6„ . 



(D7) 

(D8) 
(D9) 



apply Eq. (A7) 



Now, instead of computing the scattering matrix from 
the transfer matrix one can prove the following simple 
formulas for the scattering matrix determinant 



detS 



det S- 



+ _ 



T+T+ r p+ r p+ 

J 22^33 ± 23^32 
^22^33 ~ ^23^32 
-*ll-*44 — -*14-*41 



(D10) 
(Dll) 



Moreover Eq. (|D5| holds for e > A, while fo r e < so that det S = det S + det S . Here T?., with i,j = 1,4, 

are the elements of the transfer matrix T p . 



A is sufficient to apply the transformation Eq. (A16l 



Similarly, Eq. (D6| holds for e > |r|, while for \e\ <|T 



Defining the function 



A p (e,gL,7,0i) = e 



5 ip gi-- 



sinh 2 (6»p + ikpL) sinh 4 (^ + ik s M) + e 



— tp \qL- ' 



sinh 2 (0p — ikpL) sin 4 (k s M) 



+ 2sinh 2 (/3 + ik s M) sin 2 (k s M) [cos[(<^ T + ¥>b)/2] sinh 2 6> p + cos[27 - (tp L + tp R )/2] sin 2 (fcpi)] 



(D12) 



where, as before, k p — > ik v and p —} i0 p for |e| < |T|, 
k s —> ik s and f3 — > i/3 for e < A. The scattering matrix 
determinant is then 



Appendix E: EC free energy 



In the limit A — > the expression for A p in Eq. (D12) 
reduces to 



Al{e,qL,i,<pi) A*_{e,qL,^,ipi) 



Thus the density of states is given by 



1 d 

p(e) = - — [arg A + (e, qL, Vl ) + arg A_ (e, qL, 7, ^) 
7r ae 



(D14) 



so that 

P EC (e) 



Ap(e, gL) = sinh 2 (6'p + irfc p L) (El) 



= — J^- arctan ^coth 8 p tan sin 6 p ^j ^ 

p=± 



2L 



cosh 6L 



cos(r sinh P ) sin(r sinh 8 P ) 



r sinh 6 V 



y 

^ ± whv F cosh 2 e p - cos 2 (r sinh ( 



(E2) 
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with T = \T\L/(hv F ). Note that in the limit of large 
energy e > |T| 

4L 

lim p EC (e) = —, (E3) 



while for small energy s <C |r| 

lim p nc (e) cx . 



(E4) 



Because of Eq. (E3|, the free energy J 7 ,], defined by 



Eq. (34), diverges. However, it can be regularized by 



subtracting the constant term associated to the density 
of states at q = 0, namely p(e) = /3 E o( £ )|g=0j ancl redefin- 
ing the free energy as follows 



-t-oo 



dee[p EC {e) -/5(e)] 



(E5) 



o 



The integral in Eq. ( E5 ) can be evaluated, by intro- 



ducing a cut-off A, as follows 



1 r i r 

J-,i = — — / deep(e — hv F q/2) — — / deep(e + Hv F q/2) 

2 Jo 2 J Q 

+ edep{e) + o(^\ = 

\+hv F q/2 hv F q/2 \ 

le + J de J ep(e)+ 



2 J 

X—hv F q/2 

\+hv F q/2 



hv F q I 1 
2 



A 

hv F q/2 



de 



de p(e) + O 



A 2 



X—hv-pq/2 



(E6) 

By sending A to infinity, to evaluate the integrals above 
one needs only the limits in Eqs. (E3) and (E4) so that 



2L (hv F q) 2 ( 1 

4 V A 2 



1 



7w F 
2^ri 



This proves Eq. (43 1 



{qL? 
(E7) 



Appendix F: Bound states free energy 



Here we derive Eqs. (145]) and (46) by matching the 



solutions of the Bogoliubov-de Gennes equation at the 
interfaces between the S electrodes and the double layer 
hosting the EC. We consider < e < A so that, inside 
the superconductors, only decaying solutions are present. 
Moreover, since £ EC is much smaller than L, we consider 
only decaying solutions inside the EC. For the left inter- 
face, one needs to solve the following set of equations (in 
the linearized band approximation only the value of the 
wavefunctions must be matched and not its derivative) 



f i>Tt P (-L/2) = e^ e e^ L / 2 H mp {-L/2) 
^ Btp (-L/2) = e^^e^ 3 ^ Blp (-L/2) 
tP Blp (-L/2) = -e*P e e^-i L / 2 ^Tip(-L/2) 

[ ihip(-L/2) =e~ i ^e i ^ np (-L/2) 



(Fl) 



The system has a solution when 
expi (2p9 + 2p/3 + ip B . L - Pt,l + 2 7 - qL) = -1 (F2) 



2p/3 + 2 P 6 = ( PL -2 1 + qL+ (2n + 1)tt . 



(F3) 



Since 9 is approximately energy-independent when 
e/|r| fa A/\T\ < 1, one obtains 



arccos — p 



2 ,-parcsinl — 



and therefore 



a PL . QL . ( q£ EC 

Pl3= — - — - arcsm I 



Finally, recalling that cos j3 = e/A, 

' PL 



(F4) 



(F5) 



— = ± COS 

A 



. qL . ( q^c 



(F6) 



(F7) 



where the two solutions correspond to the two possible 
values for p. For the right interface one has the following 
set of equations 

[ ^Ttp{+L/2) = e~ l P e e^+i L / 2 ^Bt P {+L/2) 

i> Bfp (+L/2) = e^, Re -w^ Blpi+L/2 ) 

4>Bip(+L/2) = -e^P e e^+" L / 2 ^Tip(+L/2) 

{ ^ Tip (+L/2) = e~^-«e-*r^Ttp(+L/2) 

and to the related quantization condition 

expi (-2p/3 - 2p9 + 2 7 + qL + ip B ^ R - <p T ,it) = -1 

Finally, one gets 



— = ± cos 
A 



PR 

2 



qL . ( q^c 

7 h arcsm 

'2 I 2 



(F8) 
(F9) 



Eqs. (45) and (46) are obtained by observing that the 



arcsines in the cosine arguments of Eqs. (|F6[) and (F9) 



are small and amount to a simple renormalization of the 
junction length due to the finite penetration (of order 
£ EC ) of the electron wavefunctions in the EC and will be 
omitted in the following. 

An alternative derivation of the two bound states e^ 



and en can be performed by using Eq. (|D14j) to cal culate 
the density of states p BS . In the limit |T 
reduces to 



co, Eq. (D12) 
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A p (e, qL, 7, cpi) = (jv[4L-{<pt-vb)/i] + lksM ^ + e - ip [ q L-( VT - VB )/2] ^4^^) 

+ 2sin 2 (^ + ik s M) sinh 2 k s M cos[27 - (ip L + <pr)/2] 

Using that fact that 



(F10) 



sin(^ + ik s M) = e ia \J sin 2 (3 + sinh 2 k s M 

and defining 

a = arctan(cot /3tanhfc s M) = arctan(tan(7r/2 — 0) tanh k s M) , 
in the limit M — > +oo one finds that 

2Ssm[qL - ((p T - ip B )/2 + 2a] 



<p+(oc) = arg A+(e, qL, j,<p%) — arctan 



(2 + S 2 ) cos[qL - (tp T - <p B )/2 + 2a] + 2 cos[2 7 - (ip L + <p R )/2] 



(Fll) 
(F12) 

(F13) 



where 



sin 2 /3 
sinh 2 k,M 



(F14) 



is a small quantity and can be taken energy indepen- 
dent in the limit M — > oo. In Eq. ( |F13 ) we have re- 
tained term up to first order in S in the numerator and 
up to second order in the denominator. The only en- 
ergy dependent quantity left in Eq. (F13 1 is a(e). In- 



stead of calculating the density of states as a function 
of energy is much more convenient to use a as the in- 
dependent variable since what matters is the invariant 
measure (dp(a)/da)da — (dp(e)/de)de. Computing the 
derivative to the same approximation in d we obtain 



da 



where 



and 



2S (1 + cos$! cos $ 2 ) 

(cos $1 + cos <5> 2 ) 2 + 5 2 (1 + cos <I>2 cos $1) 

(F15) 



$1 = qL - ^ T - V>B + 2a (F16) 



$ 2 = 2 7 - 



(F17) 



an identical calculation for p = — gives the same result 
with a — > —a. 

This function is strongly peaked when the zero order 
term in the denominator vanishes. In the limit M — > +00 
we have a —> ir/2 — (3, so 

-cos[qL-(ip T -(p B )/2+2pl3]+cos[2'Y-(<p L +(p R )/2] = 

(F18) 

in other words 



2p(3 = 



2j-qL- (p L 



27 + qL - ip R 
so we have found the energies of two bound states 



(F19) 



± A cos (7 - 
± A cos (7 - 



qL/2 
qL/2 



Vl/2) , 



(F20) 
(F21) 



We can speak about bound states since, by expanding 
around e r or el , one can see that in the limit M — > +00 



Eq. (F15) tends to a delta function with weight it. 



We have also verified that the contribution for e > A 
of Eq. (F10) washes out for M — > +00. The only relevant 



features of the low energy spectrum are the two bound 
state peaks. 
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